
Calhoun 

iniQiuiic^iul Ar{hiv« of tilt Mil vdl Poii^roduiit School 


Calhoun: The NPS Institutional Archive 
□Space Repository 



Theses and Dissertations 


1. Thesis and Dissertation Collection, all items 


1968-09 

Digital spectral analysis. 


Dorrian, Leonard Vincent 

Monterey, California. Naval Postgraduate School 


http://hdl.handle.net/10945/40067 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 

Downloaded from NPS Archive: Calhoun 



DUDLEY 

KNOX 

LIBRARY 


htt p://w ww. n ps. e du/l ib ra ry 


Caflwuo is the Naval Postgraduate School's public access digital repository for 
research mate rials and institutiional publicatkins created by the NPS community. 
Calhoun is named for Professor of Mathematics Guy K. Caftiouo, NPS's first 
appointed — and published — schoteily author. 

Dudley Knox Library / Naval Postgraduate School 
411 Dyer Road / 1 Univefsity Circle 
Monterey, California USA 93943 






NPS ARCHIVE 
1968 

DORRIAN, L. 


DIGITAL SPECTRAL ANALYSIS 
by 

Leonard Vincent Dorrian 


UNITED STATES 

NAVAL POSTGRADUATE SCHOOL 



THESIS 

4 '- 

DIGITAL SPECTRAL ANJ&YSIS 
by 

Leonard Vincent Dorrian 


September 1968 









IIBRARY 

U S. NAVAl POSTGRADUATE SCHOOL 
MONTEREY, CALlFORK'tA 


DUDLEY KNOX I 
NAVAL POSTGfl 
MONTEREY CA 



DIGITAL SPECTRAL ANALYSIS 


by 

Leonard Vincent Dorrian 
Lieutenant, United Starves Coast Guard 
B.S., Coast Guard Academy, 1961 


Submitted in partial fulfillment of the 
requirements for the degree of 

MASTER OF SCIENCE IN ELECTRICAL ENGINEERING 

from the 

NAVAL POSTGRADUATE SCHOOL 
September, 1968 


Signature of Author 



Approved by 





Thesis Advisor 






Chairman, Department of Electrical Engin^ring 




Academic Dean 



Oo\rr\tAvq^ L. 

A program has been developed for spectral analysis using 
a medium size digital computer for computation and a general 
purpose analog computer for signal processing and system con¬ 
trol. This package is FORTRAN callable, providing the user 
with a variety of changeable system parameters in order to 
maximize the desired information required. The basic algo¬ 
rithm used is the Fast Fourier Transform. The effect of 
sampling rates, quantization and noise are investigated with 
system accuracy given. Applicational examples are shown and 
results discussed. 
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CHAPTER I 


Introduction 

Problems of Power Spectral Measurement 

The basic idea underlying all power spectral analysis is 
that a signal x(t), arbitrary to within certain limits, can 
be represented by a continuous superposition of sine waves. 
The key to an understanding of the problem is in the use and 
comprehension of the words "arbitrary within certain limits". 

It is well accepted in the Communication Engineering 
Field that, in general, the most useful relevant statistical 
properties of a signal that is considered to be stationary 
may be found in the power spectra.^ An ensemble of time 
fxonctions may be said to be stationary if any translation of 
the time origin leaves its statistical properties unchanged. 
Strictly speaking, to properly evaluate such a function one 
must measure the function over an infinitely long interval 
since the function is defined through all eternity, past and 
future. Since any measurement of the signal is necessarily 
limited in time to T seconds, we must make some assumption 
about its value outside this time. If we assiome that the 
signal is zero outside of our assxamption of T seconds, then 
the spectrum of our signal must contain infinite bandwidth. 

If we assiame that the function is periodic in the T seconds, 
then it must be composed of a discrete set of sinusoids 
either finite or infinite in number. In actual practice one 
knows he can never sample for infinite time nor can the 
physical system under measurement have infinite bandwidth. 
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Our dileirana lies in the fact that we are attempting to asso¬ 
ciate a finite spectrum with a finite time. 

Since man lives in his finite world he attempts to 
achieve finite measurements. In communications systems it 
is the changing nature of the signal that contains the intel¬ 
ligence we wish to convey. If one were to send periodic sig¬ 
nals the necessity for the use and continuance of such a sys¬ 
tem would be nonexistent. We are faced with the decision 
that the arbitrary limits imposed on our signal must be that 
it is bandlimited and of finite duration. 

If x(t) is the signal we are interested in evaluating, 
then it is customary to speak of P(f) as the power spectral 
density of the signal, where x(f) is the familiar Fourier 
integral. One deals in P(f) rather than x(f) because when 
evaluating a previously unknown signal of pseudo random com¬ 
ponents, such as phase modulated waves or admixtures of noise 

which will produce a varying x(f), x(f) may average to zero 

2 

over a T of a long duration. Thus the mathematical formu¬ 
lation is: 

P(f) = i |x(f) 1^ . 

The situation is further complicated by the fact that 
one normally assumes x(f) to be continuous, while in reality, 
the only data normally available consists of values of x(f), 
usually evenly spaced for simplicity, at a set of times, 

^1' ^2 '** ^n’ thus deals with classical time series 

rather than a continuous time function. One then in reality 
has a function that can be simulated by a finite set of 
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oscillators operating at discrete frequencies with specifiable 
amplitudes and phase relations. 

Unfortunately it is very improbable that any two samples 
of a random function provide an identical set of power levels. 
In the basic assumption of a stationary ergotic process one 
states that as T grows in length, the power determined from 
an interval, f^ to f2/ will approach the average from a large 
number of short records. Thus it is on the basis of an 
assumption of a stationary process that the credibility of 
the process of estimation of spectra rests. 

Whenever one computes a power spectrum he must have some 

3 

idea of how good the estimate is. Therein lies the final 
determination of the usefulness of any measurements. Tukey 
and Blackman^ in their famous paper state that a convenient 
description of the stability of any positive estimate may be 
found in the equivalent number of degrees of freedom such 
that: 

K = 2T * 

Through use of the Chi Square Distribution one may estimate 
the confidence interval for the measurement. This determina¬ 
tion leads to the central matter in spectral measurement, the 
uncertainty principle; that the observed fractional error of 
one's system of measurement is directly related to the mini¬ 
mum signal duration and the required spectral resolution. 

This uncertainty relationship is normally given as: 

e^T (f2-fi) = 1 

2 12 
where e is the mean square fractional error. ' 
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Need for Power Spectral Measurement 

One may ask why he may wish to obtain a complete descrip¬ 
tion of the energy spectra of a signal. In reality, spec¬ 
tral analysis is necessary, primarily from the standpoint of 
completeness of our information about the process we are ob¬ 
serving. If one desires to gain true insight, as it is pre¬ 
sumed most engineers do, one must observe both the frequency 
and the time domain. The present day need for adequate fre¬ 
quency information can be traced to the fact that it is often 
the most convenient method of summarizing data, this reason 
being directly related to the end product of the analysis. 

One concept advanced today for the use of spectral analysis 

4 

is in its use with linear systems. It turns out that most 
physical systems can be approximated by linear systems. Lin¬ 
ear Systems Analysis lends itself more readily to analysis 
in the frequency domain. Spectral Analysis thus offers a 
means of calibrating physical systems. 

It is clear that one of the most important parameters 
that must be measured in order to identify a signal or signal 
source is frequency. In any recognized adequate signal proc¬ 
essing environment spectral information is a necessity. The 
duality of the time and frequency domain demands that the 
entire content of both be explored to provide the maximum 
information content. 

Differences Between Analog and Digital Measurement 
1. Analog 

The Analog Spectrum Analyzer is an instrument de¬ 
signed to graphically present amplitude as a function 
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of frequency in a portion of the spectrum. The most 
common and best known approach incorporates a narrow- 

5 

band superheterodyne receiver and an oscilloscope. 

The receiver is electronically tuned in frequency by 
applying a sawtooth voltage to the frequency control 
element of a voltage-tuned local oscillator. The same 
sawtooth voltage is applied to the horizontal deflec¬ 
tion plates of a cathode ray oscilloscope with the 
output signal of the detector being applied to the 
vertical deflection plates. A basic diagram is shown 
below in figure 1. 
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Figure 1 
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This system has the net effect of a constant width band¬ 
pass filter with a variable center frequency. If the 
output of the detector is passed through a squaring cir¬ 
cuit, the result will be a power spectrum of the signal. 
Analog Spectral Analysis, while offering a rather mobile 
inexpensive tool for analysis, is often limited in capa¬ 
bility. This is generally a result of gradual filter 
drift and response, causing, in addition to difficult 
narrow-band filter adjustment problems, a general in¬ 
ability to perform high precision measurement. It has 
been reported that the RMS circuitry in most available 

equipment is not very accurate when repetitive impulsive 

2 

noise having a high crest factor is analyzed. 

2. Digital Techniques 

Digital Spectral Techniques., one of which serves as 
the basis for this paper, owe their present popularity 
to two events: the production of high' speed digital com¬ 
puters and the recent developments of programs for the 
computation of the Fast Fourier Transform. This allows 
a major reduction in the time*for'transforming a data 
sample from the time domain to the frequency domain. 

Only a brief description of the techniques employed will 
be given as the system will be explained in the forth¬ 
coming sections of this report. 

There are three primary digital techniques that have 
been used for spectral analysis. Each shall be considered 
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briefly. The first two techniques will be explained and 
compared with Fast Fourier Transform results. 

a. Digital Filter Techniques (without the use 
of the Fast Fourier Transform)- this technique had 
varied use prior to 1965. The method consists of 
using a group of narrow-band filters realized by 
a set of difference equations that are solved on 
a digital computer. This analysis required less 
computation than the Autocorrelation-Fourier Trans¬ 
form procedure.^ It had been reported that a 
computation time per filter had been attained of 
0.45 sec. It can be shown that for any appreciable 
analysis bandwidth and narrow spectral estimate, 
the nxmber of filters required rapidly grows. For 
example, for a 10 khz sampling rate, 1000 data 
points, and 10 hz filter bandwidth, the execution 
time would approach 250 seconds.^ 

7 

It has been shown by Kaiser and others 
that the number of weighting terms in a digital fil¬ 
ter is a direct function of the steepness of the 
filter's skirts. The amount of computation time is 
a function of the n\amber of these terms. For ex- 

7 

ample, Kaiser has shown that for a nonrecursive 
filter bandwidth of 10 hz, a sampling rate of 10 khz 
and 1% overshoot, approximately 1000 weighting terms 

g 

or "taps" are required. Stockham cites 28 as the 
point at which the Fast Fourier Transform becomes 
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the efficient method for analysis. He has shown that 

above this point, for moderate data lengths, 

savings of nearly 50 times reduction in computation 

time may be achieved. The present literature has 

8 9 10 11 

clearly demonstrated ' ' ' the superiority in cal¬ 

culation time of the Fast Fourier Transform. 

b. Semi-Hybrid Techniques: One may attempt to 
use the analog and digital computer in combination 
in an attempt to extract the best qualities of each 
for spectral analysis; bookkeeping on the digital 
side and integration on the analog side. One would 
first attempt to compute the Autocorrelation function 
R(t) and then transform to achieve the power spectral 
density function. The limiting function, in general, 
is the ability of the digital computer to output, in 
multiplex, from digital to analog form, the data 
values necessary for autocorrelation. A simplified 
diagram is shown below: 


I DlG^lTAU 



AnALOCi CoHPuTEi^ 
MULT. lUT. 


Riifi 


Semi-Hybrid Technique 
Figure 2 
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As an example, the SDS 930 is capable of 
22 micro-seconds multiplex of two channels. Thus 
for a 1000 data point set, the computation time 
would be in excess of 22 seconds for just the cal¬ 
culation of R(t) . It must be understood that the 
22 micro-seconds time quoted is only for multiplex¬ 
ing, and additional time would be required for com¬ 
puter bookkeeping operations and necessary data 
arranging. This can be compared to less than 10 
seconds for similar data using the Fast Fourier 
Transform. In reality this procedure would be one 
of substituting a digital computer for the mechan¬ 
ical operation of feeding dual data streams into 
an analog environment with the necessary delay func¬ 
tion T applied. 

c. Fast Fourier Transform: With the advent 

of this method it has been shown that it reduces 

the number of arithmetic operations from the order 
2 

of n to the order of n log 2 n. In contrast to pre¬ 
vious methods, the above discussion has shown the 
savings in computational time that may be achieved. 
This improvement has led to basic conceptual re¬ 
visions in computing methpds^ for spectral analysis. 

The first step in spectral analysis,-is "to dig¬ 
italize the analog•signal by- sampling at a rate that 
is twice the highest frequency present. Generally 
one is then left with the choice or proceeding in 
the time domain or in the frequency domain as shown 
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in figure 3. The fastest route to the Autocovariance 
Function is now usually through the Fourier Trans¬ 
form and Inverse, yielding the power spectrum as a 
by-product. 

This transform method will be covered in detail 
in a later section of this report. It is readily 
seen from the previous comparisons that this method 
offers the most rapid method of analysis presently 
available. It has been chosen as the basic method 
for analysis in this report in order to demonstrate 
the capabilities, limitations and possible exten¬ 
sive use of this new method. 

Concept of a Callable Frequency Analysis Routine 

The USNPS computing facilities include an "open shop" 
computer systems laboratory operated by the Department of Elec¬ 
trical Engineering. This laboratory includes a medium size 
digital computer and a general purpose electronic hybrid/ 
analog computer. Soon after the installation and acceptance 
of this new facility in late 1967 it became apparent that the 
use of the facility would be varied and available for inter¬ 
disciplinary use. The need for callable siibroutines to per¬ 
form desired functions is obvious to all who have any famil¬ 
iarity with general computer facilities. 

This writer became interested in modern speech research 
as a possible study area. One important mode of that field 
today is concerned with spectral analysis. It became appar¬ 
ent that a rapid accurate means of digital spectral analysis 
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was a necessity. Due to the nature and use of the laboratory 
it was decided to attempt to formulate a spectral analysis 
procedure that was FORTRAN callable. It was also to be made 
available for interdisciplinary use with the minimum of pre¬ 
liminary investigation and study required by a potential user. 

A spectral analysis procedure was developed and is the 
basis of this report. It is felt that this package will offer 
a valuable means of investigation to any potential user inter¬ 
ested in spectral analysis and will lend itself to use as a 
basis for further study and application. A detailed descrip¬ 
tion, theory of operation, results and conclusions obtained 
will be given in the'following sections of this report. 
Description of Facilities Used 

'The Computer Laboratory of the Electrical Engineering 

Department contains a mediiam size digital computer, Scientific 

Data System-930, and a general purpose hybrid/analog computer 

Comcor Inc. 5000. Figure 4 shows the present and planned 

12 

facilities of the laboratory. 

The SDS 930 computer is a high speed, general purpose 
digital computer with the following characteristics: 24-bit 
word plus parity; binary arithmetic; one bank of cell memory 
16K; 1.75 micro-second memory cycle time; fixed point add 3.5 
micro-‘second and 7.0 micro-second multiply; floating point 
add 14 micro-second; multiply 123 micro-second and an exten¬ 
sive software capability. The versatility of the system is 
exemplified by its input-output capabilities. Input/output 
devices include teletype paper tape, punched cards and digital 
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Figure ~4 
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magnetic tape units, a continuous function generator and con¬ 
tinuous signal magnetic tape inputs through the analog-digital 
interface from the analog computer. 

The Cl 5000 is a large scale, general purpose electronic 
analog/hybrid computing system. The system has a highly de¬ 
veloped analog/hybrid control and functional capability. 

In the package developed, the time varying signal of 
interest is input into one of the Analog-to-Digital conversion 
trunklines of the Cl 5000. This converter has a 14 bit con¬ 
version capability with an aperture time necessary to provide 
a maximum sampling rate of approximately 50 khz. A priority 
interrupt capability is used to determine the actual sampling 
rate by using the master clock in the logic section of the 
Cl 5000. Standard frequencies of 100 khz, 10 khz, 1 khz, 

100 hz, 10 hz, 1 hz and 1/4 hz are available. In addition, 
through the use of the logic board, preset values of frequency 
may be achieved by chaining logic sections to generate arbi¬ 
trary frequencies. In this manner th\ambwheel counters may be 

set and a frequency determined by dividing the values of the 

13 

preset covinters into the standard clock frequency selected. 
Positive program control is maintained through the use of the 
sense line capability of the SDS 930. The input signal and a 
bias threshold level are applied to a comparitor in the analog 
side of the Cl 5000. When the output of the comparitor 
achieves a non-negative value, a pulse is sent to a delay flip 
flop on the logic board. This sets a positive value on the 
sense line. The digital analysis program will wait in an idle 
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loop until the sense line indicates a positive logic. In this 
manner the user has the capability of placing the program in 
run and commencing data input and analysis at his discretion. 

The magnetic tape xinits were used as a temporary storage 
device and a ready interface between the FORTRAN II format 
and machine language programming necessary for high speed 
sampling. This necessary relationship will be covered in the 
next section of this report. ‘ Digital Analysis was performed 
on the SDS 930. Processed spectral information was outputed 
on the line printer in graphical format for ease of analysis 
and high resolution capability. A sense switch on the SDS- 
930 was used to allow the user to have the data values of the 
spectrum printed if it was so desired. 
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CHAPTER II 


DESCRIPTION OF PACKAGE DE^/ELOPED 
Basic Program Context 

The concept of versatility provided the main criteria 
for the development of this program. An attempt was made in 
the calling argxjment to provide the potential user with the 
necessary parameter variables so that he would have the abil¬ 
ity to change the method of analysis and observe the effects. 
It was felt that this would allow the program to be used for 
more potential applications than one that was restricted to 
one mode of operation. 

A listable output of the program as well as the META¬ 
SYMBOL subprogram is listed in Appendix A. The program calls 
a subprogram HARM for the purpose of computing the Fast 
Fourier Transform. 

The subroutine DATA is called in the following manner: 

CALL DATA (FREQ, TIME, MBITS, Ml, SMOU, Fl, F2) 

Description of parameters: 

FREQ - This determines the sampling rate that will be 
used. It is the clock pulse that is generated in 
the logic side of the Cl 5000 and inputed into 
trunkline T210 for use in the priority interrupt 
system. 

TIME - This specifies the length of time that the signal 
will be sampled. Due to core size limitations it 
is necessary to maintain a limit of 1024 on the 
empirical relationship: 

FREQ*TIME*NBITS/24 < 1024 
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NBITS - This specifies the amount of bit capability that 
will be used from the output of our 14 bit convert¬ 
er. Artificially quantizing by not using the full 
resolution capability of the converter allows the 
program the provision of increasing the number of 
samples taken. 

Ml - This specifies the data size of the Fast Fourier 
Transform that will be used. The use and capabil¬ 
ity for different size data inputs to the nrogram 
will be explained in the next section of this re¬ 
port . 

SMOU - This specifies whether or not a smoothing func¬ 
tion would be used on the final Fourier coeffi¬ 
cients generated; SMOU = 0, no smoothing function 
used, SMOU = 1, a smoothing function is used as 
follows: 

'■k = T '-^'k-l+^k-^k+l' 

This procedure is prevalent in the literature^^' 
Singleton and Poulter explained this as an addition¬ 
al form of digital smoothing, performed with the as¬ 
sumption that the power spectral density function 
is continuous and thus used as a form of averaging 
across the frequency band to gain statistical sta¬ 
bility. The ability to choose whether or not to 
use this smoothing was felt to be an important tool 
in the hands of the user. This effect and its use 
were noted and will be covered in a later section. 
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Fl, F2 - These parameters are included as a means of ex¬ 
panding the spectral plot to offer the analyst the 
ability to display an area of particular interest 
with a very high degree of resolution in addition 
to an overview of the complete bandwidth. 

Use of Machine Language-Fortran Mixture of Subprogram 

The use of a META-SYMBOL (SDS Machine Language) sub¬ 
program in this package was based on the requirement to 
achieve the maximum sampling rate possible. Under the SDS 
Real-Time Monitor system delivered with the 930 system, the 
F-IV system was 'capable of a maximum sampling rate of 2 khz. 

It was felt that this sampling rate was well below the accept¬ 
able lower limit deemed necessary for adequate analysis in 
the audio range. 

A machine language sxibprogram was developed with a wide 
range of capability to act as an important part of the anal¬ 
ysis package. This subprogram provides the desired versa¬ 
tility required by allowing several of the parameters of the 
sampling process to be callable. A maximum sampling rate of 
12.5 khz was achieved prior to the onset of aliasing, the 
folding back of frequency information that occurs when the 
Nyquist criteria of sampling at twice the highest frequency 
present is not met. 

The subprogram Input is listed in Appendix A as part of 
the main subroutine DATA. The subprogram is called in the 
following manner; call INPUT (INPl, INP2, LENGTH, NUMRCDS, 

N). Description of parameters; 
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INPl - This is an argioment necessary for bookkeeping 

operations and the generation of certain shifting 
and masking operations. For ease of calculation 
and simplicity, the value of these parameters are 
calculated in the main stibroutine DATA. 

INP2 - This has the same function as INPl. 

LENGTH - This is the total ntimber of samples to be taken 
less one. The value is reduced by one because of 
its particular use and relationship to a machine- 
language cycling operation and test, SKR. (Appen¬ 
dix B gives a listing of the machine-language 
operations and their meaning in reference to the 
basic register manipulations of the SDS 930). 

NUMRCDS - This is a necessary argument used in the book¬ 
keeping operation required to write data on magnetic 
tape in a format compatible with Fortran II opera¬ 
tions. The data is arranged in blocks of a 

specified length with check sums to insure data 
accuracy and completeness. A section of the machine 
language sxibprogram was devoted to acting as a 
servicing agent for the sampling program. This 
acted in a form that one might call an artificial 
buffer. 

N - This is the same parameter as was listed in the 
subroutine DATA as NBITS. 
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Discussion of Program Flowqraphs 

Flowcharts for subroutine DATA and for subprogram INPUT 
are shown in Figures 5 and 6. A flowchart is provided for 
subprogram INPUT to offer insight into the different processes 
that were required for proper data handling in the sampling 
procedure. 

It is felt that the flowchart for DATA is self explan¬ 
atory; however there are a few portions of INPUT that are 
worthy of note.' The program constants calculated in step 2 
are generally required for the various counters and masks 
necessary for data manipulation. In addition, subroutines 
HARM and PLOTP require various entering arguments. These con¬ 
stants are all interrelated and are determined once the nvim- 
ber of data points in the Fast Fourier Transform (HARM) are 
chosen. They are also used to determine the proper values of 
the masks and counters to be used with the non-standard-length 
words for data storage. 

The type of data manipulation is determined by the bit 
size chosen. This is directly related to the packing process 
shown in the flow chart. For example if a bit size of 12 was 
decided on we would then pack two samples into one storage 
cell as shown below. 


Data word 1 

111 

Ill 

Tor 

101 

110 

000 

obo 

000 

Binary 


7 . 

7 

5 

5 

5 

0 

0 

0 

Octal 


Data word 1 

Quantized 7755 


(Note: using only the 12 most significant bits) 
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Data word 2 

Oil 

100 

100 

111 

010 

000 

000 

000 . 

Binary 


3 

4 

4 

7 

2 

0 

0 

0 

Octal 


Data word 2 

Quantized 34 4 7 

(Note: using only the 12 most significant bits) 

Resultant 
Packed storage 

Cell 77553447 

In order to accomplish this, masks must blank the unwanted 
bits. Special shifting counters are used to arrange data 
words in the proper order when placed in the final storage 
location. In a system where core storage is at a premium, 
if one can allow for quantizing to exist in certain applica¬ 
tions , a valuable saving of storage areas can be effected. 

For example, of one is satisfied with a bit size of only 6, 
it is then possible to have an effective core area four times 
as large as that used for a 24 bit word. 

Once the data has been acquired it is then necessary to 

write this information on magnetic tape in standard format 

X 5 X 6 

compatible with FORTRAN II. ' Data manipulation is thus 
performed in the machine-language subprogram, the only possible 
solution to acquire the needed flexibility. The unpacking 
process is the reverse procedure to that followed in packing 
the storage cell. This data is in fixed point format. After 
input into the main subroutine it is necessary to convert this 
information into floating point to be compatible with the Fast 
Fourier Transform subroutine HARM. 
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Timing Restrictions 

The cycle time of a computer is that time reference unit 
normally used as the basis for judging the relevant speed of 
operation of all execution cycles of the machine. Once the 
design of a computer is completed and the cycle time is fixed, 
there is then an absolute lower limit placed on the speed of 
operation of basic^mathematical manipulations. It is the 
ingenuity* of the programmer that determines how close one can 
approach this lower limit by careful and precise use of ma¬ 
chine instructions. One can judge, by the machine cycle time 
for example, the ability of a particular computing system to 
perform rapid analog to digital conversion and the subsequent 
storage of this information. 

A controlling parameter on the procedure developed was 
the maximum rate at which data samples of analog to digital 
conversion could be taken and stored in their final cell 
locations. It was determined that the maximum rate possible 
was 12.5 khz when one retained the capability for bit dis¬ 
crimination and svibsequent cell packing. 

An understanding of the Priority Interrupt System used 
in the SDS 930 computer is vital if one is to have insight 
into the limitations imposed on the maximum sampling rate. 

The SDS 930 Interrupt Control System "provides added program 
control of input/output and compute operations and allows 
immediate recognition of special external conditions." An 
interrupt is a signal external to the main computer frame 
which causes the computer to transfer program control to one 
of a selected set of memory locations. Both the interrupts 
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used in the sampling program were generated external to the 
digital computer; interrupt 0200 by the Analog—to-Digital Con¬ 
verter and interrupt 0210 by the clock in the logic side of 
the COMCOR 5000. When an interrupt occurs (a pulse on the 
interrupt wire) the computer executes the instruction in the 
memory corresponding to the interrupt wire. An instruction 
transferring control to a specified servicing routine was 
placed in cell locations 0200 and 0211. 

As the name implies, the interrupts have a priority sys¬ 
tem. Interrupts are arranged in order of precedence by numer¬ 
ical hierarchy. The highest priority is always served even 
though another interrupt may be in the waiting state, that is, 
when the interrupt signal has been sent but the computer is 
servicing an interrupt of higher priority. Interrupts are in 
an inactive state with reference to the main computer when 
they are disabled. It is necessary to enable interrupts (a 
machine instruction) prior to any operation requiring inter¬ 
rupt control. It is also necessary to service any interrupt 
in the waiting state prior to disabling interrupts in order 
to prevent the control of the program from being inadver¬ 
tently lost to an interrupt routine from a previous portion 
of the program. Inattention to this concept of control will 
result in some surprising results and inherent frustrations. 

Listed below is the central portion of the sampling 
routine: 
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Instruction 


Niomber of cycle times 



EOM 

034000 

1 


POT 

cw 

3 


EOM 

033004 

1 


EIR 


1 

CYC 
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1 
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BUR 

$-2 

1 

AD IN 

PZE 

0 

1 


EOM 

034000 

1 


POT 

CW 

3+ 


BRU 

*ADIN 

2 

SAMPLE 

PZE 


0 


DIR 


1 


LDA 

LOC 

2 


ETR ■ 

MASKl 

2 


MRG 

*COUNT 

3 


BRU 

$+2 

1 


BRU 

$+4 

1 


LRSH 

0,2 

7 


STA 

*COUNT 

4 

BRU 

BRU 

$+5 

1 


STA 

*COUNT 

4 


MIN 

COUNT 

3 


LDA 

INP2SM 

3 


STA 

PARAM2 

3 


SKR 

SMLENG 

3 


BRU 

*BACK 

2 


BRU 

*$+l 

2 
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Since we are interested in the maximum sampling rate ob¬ 
tainable, the nimber of cycle times required for the service¬ 
able routines are calculated for the worst case in a particular 
mode, that is, the one path depending on the option that takes 
the longest. 


Basic Loop 

3 

ADIN 

6x2 

SAMPLE 

30 

TOTAL 

= 45 cycles x 1.75 micro-seconds 

T 

= 78.75 micro-seconds 


fm = 1/T = 12.6 khz 

Note that the A/D routine ADIN was multiplied by a factor of 
2 because one could receive an interrupt 0200 (ADIN), and 
0210 (SAMPLE), simultaneously and ADIN would thus be served 
first. The cycle could conceivably be ADIN, SAMPLE and then 
ADIN prior to the next sample being taken. It is to be noted 
that this theoretical result agrees closely with the observed 
result of 12.5 khz. Aboye that frequency of sampling the 
program begins to miss samples. ' . 

Type Results Observed 

Figure 7 is an example of the type of results observed. 
PLOTP, the plotting subroutine used, gives the x axis coor¬ 
dinates in linear form and the y coordinates in logarithmic 
form. This provides the maximum range of presentation pos¬ 
sible. The nonlinearity observed in the plot usually will 
tend to detract very little from the informational content of 
our presentation, since the area of greatest nonlinearity. 
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the lower amplitude level, contains noise information. It is 
normally noise that we wish to suppress in order to emphasize 
our signal information. The plot shown is of a sinusoid at 
3.001 khz. 
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CHAPTER III 
THEORY OF OPERATION 
Analog to Digital Conversion 

The conversion of signals from analog to digital form is 
generally a two-step process. First, the analog signal of 
figure 8 (a) is sampled in time, figure .8 (b) , the samples are 
quantized in amplitude as shown in figure 8 (c). 


Cd) 

(C) 



/)a<ALOG, SiG)MAU 

Sampled 



Sampled 

QuAMTlltS) S^fJAL. 


Figure 8. The Sampling and Quantizing Process 
18 

Bennet , in his often quoted work, characterized samp¬ 
ling as "quantizing of time". Sampling is the process of 
attempting to specify a continuous signal by assigning dis¬ 
crete values of amplitude for specified positions of time in 

such a manner as to completely define the signal. Nyquist 

19 

and later Shannon set forth the basic rules that govern the 
sampling theorem: 

1) Every signal of finite energy and bandwidth W hz may 
be completely recovered from a knowledge of its samples taken 
at the rate of 2W per second. 
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2) Such a function, bandlirtiited and of finite energy 
may thus be said to comprise a region in an infinite dimen¬ 
sional "function space" characterized as being limited to a 
time T in a space of 2 TW. 

To the unwary, the Nyquist rate and Shannon's concept of 
Signal Space may seem a completely definitive answer for the 
sampling question. However one must realize that this was 
offered as rather an absolute bound: signals exist neither 
in bandlimited form and for all eternity nor in a finite time 
and infinite bandwidth. It is often said that sampling at 
higher than the Nyquist rate offers no further information 
about the. signal. This would be true, were it not for the 
presence of noise and also for conversional errors often mis- 
characterized as noise. One often hears of quantizing errors 
and sampling errors being designated as analogous to various 
forms of noise. 

Conversion errors are unavoidable since we are attempting 
to translate a continuous world into a digital form with con¬ 
tinuous input, continuous output and continuous relations be- 

20 

tween variables. The sources of these errors are: 

1) Sampling errors: Time is sampled and the input sig¬ 
nal is considered at sampling instants only. 

2) Quantization errors: Every variable has to be repre¬ 
sented by a finite number of values. This results in the 
rounding off of the continuous value of an input variable 
and the generating of round-off errors whenever a com¬ 
putation is performed. 
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In each case a quantity known to great precision is expressed 
with considerably less precision. 

It is the finite limitation on time in the sampling proc¬ 
ess that leads to sampling errors in time. This is often 
called aperture error or frame error. It is that amount of 

time required by our physical sampler to make the decision 

20 

as to which level our signal is at. Katzenelson has char¬ 
acterized this as a sampling efficiency problem defined as: 

, -2aT 

- 

^s 2aT 

where "a" is a constant of the probability density function 
of the input and T is the length of time required to sample. 

It is seen that as T ->• 0, approaches 1. But it must also 
be recognized that physically, T can only approach zero and 
will never result in 100% efficiency. 

In the quantizing process we nomally divide the maxi¬ 
mum amplitude range of the signal (Xq,X^) into L non¬ 
overlapping subintervals such as {X 2 i^ 2 ^ • ^ finite L 

some error must be incurred. The obvious tendency would be 
to attempt to make the grain (fineness of the interval) as 
small as possible. The attempt to decrease grain size, how¬ 
ever, is immediately met with increased cost, as the number 
of intervals required increases. It is well understood that 
when the grain size decreases, the precision of the measure¬ 
ment increases. However an important effect to be noted is 
that the confidence of out measurement proportionally decreases 
This is due to the inescapable fact that noise always pre¬ 
vents 100 percent confidence. 
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The sampling rate, that is the nvunber of samples we wish 

to take per cycle, has been shown to have an effect on the 

21 

confidence of our measurement. Initially one must insure 
that the Nyquist criteria are met. If the signal of interest 
X(t) is available for uniformly spaced intervals of t. At, and 
we have assiamed that the signal is bandlimited, the notion of 
"aliasing", the so called "strobing effect", enters if there 
. are, in fact, contents of the signal above the band limit 
point f = 1/2 At. All frequencies above f are indistinguish¬ 
able from the bandlimited spectrum of 0 - f. Blackman and 
Tukey have called this the "aliased spectrum". Figure 9 shows 
that when one fails to observe this interrelationship between 
the bandlimited condition auid the Nyquist rate, one slows down 
and folds back frequency content, leading to erroneous con¬ 
clusions. 




Figure 9. Aliasing Effect 
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Sine Waves whose frequencies differ by multiples of 1/At can¬ 
not be distinguished and are folded back into the band 
(0 - 1/At). 

The sampling rate has been found to have a direct effect 

on the quality of fidelity of the reconstruction of the sig- 

21 

nal from the discrete samples. Melton has derived an expres 
sion to indicate the effect of the sampling rate upon the 
error that may be increased when one samples in the worst 
possible time condition of the cycle. His error expression: 

E = l/[cosII/n-l] 

where n is the number of samples per cycle, indicates that 
for a certain grain of quantization one may conceivably in¬ 
cur infinite error as the expression "blows up" at n = 2. 

This expression does indicate that using more than 8 samples 
per cycle yields little improvement in fidelity. 

It has been stated that with the presence of noise, one 
can never be assured of achieving 100 percent confidence. A 
higher rate of sampling does not reduce the effect of noise. 
The critical relationship between noise and sampling is that 
the sample points must be independent with respect to noise; 
that is, the noise voltage at one sample point must not in- 
fluence the noise voltage at. the next. This is just one way 
of stating that the Nyquist criteria also apply to the noise 
present. One must sample at a rate at least twice as high as 
the highest component of noise present. 

The number of binary digits required to reconstruct an 
analog signal depends on the noise present and becomes crit¬ 
ical at the lowest signal levels, where the problem becomes 
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one of detection of the signal. If the noise is considered 

Guassian in form, a good rule of th\imb, based on Chebyshev's 

Inequality, states that for a confidence level of approximatelv 

75% the quantization grain should be set at approximately 3 

21 

times that of the RMS value of the expected noise. 

Thus by use of the bit size of the analog-to-digital 

converter one may in fact control the quantization size. 

18 

Bennet has shown that as we increase the bit size above 8 
digits the distortion spectrum is relatively flat and concen¬ 
trated in the low-frequency area. Quantizing error due to 
distortion without the presence of noise is a major factor to 
be considered when one is attempting to minimize the number 
of digits used. 

Thus in choosing our grain size and the number of bits 
used, the dual factors of inherent error noise and the rela¬ 
tionship of the grain size to imposed noise must be considered. 

When a converter has a capability for a set number of 
digits and one chooses to use only a portion of this number, 
the effect is then to increase the grain size as well as 
decrease the number of bits. One then is forced to consider 
both noise requirements, induced and imposed noise. 

Results of the Analog-to-Digital Conversion Program 

The results obtained by the sampling process were all 
performed on sinusoidal signals with a constant number of 
points (512) in the Fast Fourier Transform subroutine and 
without the use of a smoothing function on the final Fourier 
coefficients. The sinusoidal input signals were chosen as 
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a basis for the analysis due to the well known Fourier Trans 
form; 

3 ^ [cosojQt] = j[6(f-fQ) + 6(f+fQ)] 

where 6 is the Dirac Delta Function. It was felt that this 
would give the greatest precision in measurement as well as 
the maximum signal-to-noise ratio attainable and offer the 
best use as a standard. 

Full bit-size capability was used in all cases except 
when analyzing the effect of altering the bit size. This was 
done to isolate, as much as possible, the interrelationship 
that results from parameter changes in the sampling program 
and that of the data analysis program. One may consider the 
sampling program'as ideally separate. It is merely a method 
of generating data for use in the analysis program. The anal 
ysis program and subsequent output plot does however offer an 
excellent vehicle for the study of changing parameters and 
resultant effects on the quality of sampling attained, A 
typical plot is shown in Figure 10. 

1) Aliasing Effect 

Figure 11 is the resultant frequency spectrum plot 
of a sinusoid of 6.5 khz sampled at 10 khz, measured 
under identical conditions to the signal shown in figure 
10. The well known "folding back" effect is obvious. 
With a Nyquist Rate of 10 khz our bandlimited region is 
0-5 khz. At 6.5 khz we are 1.5 khz above the region 
and our "aliased" spectrum indicates the erroneous 
position of 1.5 khz below the band limit. Figure 12 
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shows this effect in a plot of true frequency versus sam¬ 
pled frequency for sinusoids taken from 1 khz to 9.5 khz 
at a 10 khz sampling rate. In each case shown a Hewlett- 
Packard 5245L digital counter was used as the reference 
standard. Its tested accuracy is one bit in the last 
decimal place shown. 



Figure 12 

2) Effect of Varying Bit Size 

Figures 13 through 18 indicate the effect of vary¬ 
ing the bit size of the sample from 1 bit to 24 bits, on 
a 1 khz sinusoid scimpled at 10 khz. It must be noted 
that the Analog-to-Digital converter used is a 14-bit 
converter so that any bit level above 14 does not add 
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further information content to our sample. In the plots 
shown a maximum input of 90 volts was used. 

This was done in an attempt to use the full capability of the 
converter and to input the maximum signal-to-noise ratio pos¬ 
sible. The value of each bit in the converter is shown below: 


Position 

Value 

0 

sign bit 

1 

> 

O 

in 

2 

25 V 

3 

12.5 V 

4 

6.25 V 

5 

3.125 V 

6 

1.562 V 

7 

0.781 V 

8 

0.390 V 

9 

0.195 V 

10 

0.097 V 

11 

0.048 V 

12 

0.024 V 

13 

0.012 V 


Table 1 


One can thus see that the precision of the measurement is spec¬ 
ified by the bit level chosen. Specifying a bit level of 6, 
which would be using all the bits up through bit 5, determines 
a precision of ±3.125 v. All signals above or below the value 
of 96.875 V would have the same value. We have thus chosen a 
grain size of -3.125 v and the arguments set forth in the first 
part of this chapter apply. 
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For a sinusoidal input of sufficient magnitude to 
full capability of the converter very little effect can 
noted on either the frequency or the amplitude. In the cas 
of a bit size of 1, it can be seen that we are in actuality 
measuring only the zero crossings of our signal and that we 
introduce even harmonics. This can be easily shown through 
the standard Fourier Analysis of a square wave as given below. 








_^ 

L 


ZE ^ 


-('u^ = 1 2 j. SIN i][y: 
IT <1 L 


By our quantizing scheme we thus are entering frequencies into 
the analysis which may violate the Nyquist rate. 

This effect diminishes, however, above a bit level of 3. 
The effect that does predominate is one of an artificial bias 
inserted by our quantization scheme. Due to the binary con¬ 
figuration of the words of the computer the use of two's com¬ 
pliment arithmetic imposes a bias on the system. Note that 
for a bit level of 3 a signal of ±12.5 volts magnitude would 
be indicated by a slight negative bias. When the sign bit is 
lit negative the binary representation must begin with a "1". 
This is converted into octal and we are faced with a very 
large negative number, 40000000, which is input as a data 
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sample into the analysis routine, it can be noted, however, 
that this effect decreases sharply as the bit level increases. 

The signal-to-noise ratio, however, does decrease slight¬ 
ly as the bit size decreases. This is due to the fact that 
as we increase our bit size we decrease our grain size. This 
allows the sample value to indicate more precisely small ran¬ 
dom fluctuations. Thus the presence of noise is indicated 
with an increase in noise level. It is surprising to realize 
that the course grain of quantization actually acts as a 
rough digital filter. If however, the signal has not used 
the full capability of the converter, an increase in bit size 
would not have noticeably improved the signal-to-noise ratio. 
It is actually a result of the cross effect of quantization 
noise and system noise. For quantization noise, as the grain 

size decreases the mean square of the noise decreases as: 

2 2 
n = LV12 

where L is the grain size. 

3) Filter Effect on Non-Bandlimited Signals 

In order to insure that our signals were bandlimited 
a Krohn-Hite 330-A active filter was used prior to the 
sampling. Other means were available, such as using 
active operational amplifiers on the analog computer, 
but the Krohn-Hite was chosen because of its simplicity 
of operation and the temporary lack of equipment such as 
circuit boards in the laboratory facility. 

Figure 18 is the frequency spectrum of a 30 v square 
wave at 1.025 khz. Its frequency content is in agreement 
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with the Fourier Analysis of a square wave given earlier 
in this chapter. One may notice that if the sampling 
rate used had been 4 khz instead of 10 khz considerable 
misinformation would have resulted due to "aliasing" of 
frequency content at 3 khz and 5 khz. The result clearly 
indicates the effect of the filter and demonstrates the 
necessity of bandlimiting if one is to give true fre¬ 
quency content in the interval chosen. See Figure 19. 

4) Effect of Exceeding Maximum Sampling Rate 

As noted previously in this report the maximimi sam¬ 
pling rate attainable was approximately 12.5 khz. Fig¬ 
ures 20 through 22 show the results in a plot of the 
frequency spectrum of sampling a 1.0 khz sinusoid at 12.5 kh^ 
14.3 khz, 16.668 khz and 20.0 khz. This was done to 
ascertain the effect that is achieved when the maximum 
sampling rate is exceeded. 

It is seen that in general the effect is one of 
almost an inverse of "aliasing". It actually makes the 
signal appear at a higher frequency than is the case. ' 

An explanation is more easily achieved with the aid of 
Figure 24. For the purpose of clarity and ease of ex¬ 
planation it will be assumed that the maximum sampling 
rate attainable was 10.0 khz. Figure 24 shows a 1 khz 
sinusoid sampled at 10.0 khz.'' When we attempt to sample 
at a rate higher than the maximum rate we continue in 
' real time to sample at 10 khz. The computer is not under 
real-time control and assumes that it is receiving sam¬ 
ples at a rate of 20 khz. This is due to the fact that 
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the method for control of the length of time sampled was 
one of using a counter. Each time the sample routine 
was entered it was decreased by one. This original 
counter value was the number of samples we would take. 

It was equal to the length of the sampling time multi¬ 
plied by the sampling rate assumed. Both sampling 
rates input the same data points; however, the 20 khz 
program assumes the samples are at one half the time 
spacing, thus the erroneous conclusion is reached that 
we have twice the value of the frequency present. This 
is shown in Figure 24. 

Fast Fourier Transform 
1) General 

The Fast Fourier Transform (FFT) is a method for 

efficiently computing the Discrete Fourier Transform 

(DFT) of a time series (discrete data samples). This 

transform method was brought into its present state of 

exceptionally popular use and acclaim by Cooley and 

Tukey in 1965. Most historians credit the early work 

in the modern Fast Fourier Transform to Runge who in 

1903 described a computational technique for 12 and 24 

point Fourier Transforms.Without the aid of high 

speed computational services and rapid accurate devices 

for attaining discrete data, interest in this approach 

remained purely academic. Some, work was done on this 

23 24 

approach and various similarvO.nes ' in the late 1930 s 
and early 1940's but it was not until the Cooley and 
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Tukey paper that a renewed interest in the Fourier Trans¬ 
form as an analysis tool was felt. The main reason for 
this high interest was that standard computer programs 
in use were using up extremely long and thus inherently 
expensive computer time to calculate the Fourier Trans¬ 
form of discrete data samples. The amount of computa¬ 
tional time required for standard programs was propor- 
. 2 

tional to N , where N is the number of data points, 

while the Cooley and Tukey method reduced this amount to 

N log 2 N computations. This new-found efficiency results 

in savings of computational time such that programs on 

an IBM 7094 of 8192 samples may now be calculated in 

under 5 seconds while under previous methods times in 

2 ^ 

excess of 30 minutes were not uncommon. 

The Fast Fourier Transform is thus a highly effi¬ 
cient method of computing the Discrete Fourier Trans¬ 
form of a time series. ' The use of the Discrete Fourier 
Transform is a necessity if digital spectral analysis 
techniques are to be used for the analyzing of contin¬ 
uous waveforms. This transform has properties entirely 
analogous to the Fourier Transform of a continuous wave¬ 
form if the basic criteria for sampling as given by 

25 

Nyquist and Shannon are met. It has been shown that 
the Discrete Fourier Transform of a sequence of Nyquist 
samples is directly related to the Fourier Transform of 
a continuous waveform where the Discrete Fourier Trans— 
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form is defined as 



'! 

n-l - 

A = I X, EXP [-2IIjrk/N] r = 

^ k=0 ^ 

where is the coefficient of the DFT and X^^ denotes 

the sample of the time series which consists of N 

samples. 

The transform properties of the DFT offer an impor¬ 
tant tool in spectral analysis.One may convolve two 
time series of sampled continuous inputs when faced 
with the product of the two discrete transforms. The 
shifting property is also applicable, with slight mod¬ 
ification, from the continuous transform and the DFT of 

the svici of two functions is in fact the sum of the DFT 

25 

of the two functions. 

As mentioned previously it is the ability of the 
Fast Fourier Transform for rapid calculations that leads 
to its widespread use and application. It not only re¬ 
duces the computation time but also reduces the roundoff 
errors associated with the calculation by its drastic 
reduction in total calculations. It is basically a 
"computational technique of sequentially combining the 
DFT of individual data samples such that the occurrence 
times of these data samples are taken into account se¬ 
quentially and applied to the DFT's of progressively 
larger mutually exclusive subgroups of data samples, 
which are combined to ultimately produce the DFT of the 
complete series of data samples. 
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^ - Package Obtained for This Program 

The s\ibroutine used to calculate the Discrete Fourier 
Transform used in this program is the standard version 
of the Cooley and Tukey Algorithm as set forth in "An 
Algorithm for the Machine Calculation of Complex Four- 
ler Series . This version was modified for the SDS 
930 FORTRAN-II by Professor R. L. Johnson, U.S. Naval 
Postgraduate School, with only minor modifications to 
make the program compatible with the capabilities of 
the FORTRAN package used. An outline of the subroutine 
used and the basic algorithm is contained in Appendix 
C. 


The speed of computation observed was approximately 
4 seconds for the 512 point version and 10 seconds for 
the 1024 point version. The speed of operation was 
limited by the lack of floating-point hardware. 

It must be understood that the Discrete Fourier 
Transform gives frequency content at discrete points in 
the spectrum only. These points are related to the sam¬ 
pling rate, the number of points and the frequency band 
over which they are taken, as shown below: 

(Sample Time) DT = 1/FREQ 

where FREQ is the sampling rate 

(Band Limit)' BW = 1/DT*2 

(Frequency of DEL = 1/DT*NPTS 

Measurement) 

where NPTS is two times the 
number of data points 

Figure 25 shown below demonstrates these relationships: 
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Figure 25 

Thus it can be seen that for an input of 512 data 
points we determine 513 spectral estimates. Half of 
these are negative estimates and one is the zero or the 
D.C. component. These points are spaced at intervals 
of (^)fs/256 and represent the square root of energy 
contained in the band of width DEL centered on the 
point. 

In order to narrow the band "DEL" of our measure¬ 
ment it is necessary either to reduce our sampling fre¬ 
quency or to increase the number of points used. 

Usually one wishes to maintain as high a rate of sampling 
as possible. The normal path is to increase the number 
of points used. As an example^DEL, the spectral estimate 
for a sampling rate of 10 khz and 512 data points is 
approximately 19.53 hz. 
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3) Use of Smoothing Functions 


Pi^®valent in the literature is the necassity of 
applying certain "windowing effects" to prevent the 
interaction of spectral estimates. Blackman and Tukey^ 
introduced this idea in 1948 with the concept of form¬ 
ing moving linear combinations with unit total weight. 
This smoothed function resulted in the weighted aver¬ 
ages of nearby values of the original function being 
summed in order to produce a new function. The main 
reason for using the various windows was to increase 

the stability over the spectral estimate and to pre- 

2 8 

vent leakage. Brigham has defined leakage as "that 
phenomenon by which the value of the spectral estimate, 
which was intended to estimate a weighted average of 
the spectrum over a specific frequency interval, is 
affected by components of more or less remote frequen¬ 
cies outside the interval." Many basic windows have 
been used such as: 

1) Bartlett Window; 

T (sinnfT /nfT^)^ 
m m' m 

where T is time between samples, 
m 

2) Hanning Window; 

0.5 + 0-25[Aj^_3_+Aj^+i] 

\^here A, is the Fourier coefficient. 

3) Hamming Window; 



4) Many variations, such as a basic Bartlett Window 
multiplied by constants and the sin x/x function raised 
to higher powers to steepen the slopes cf the window 
in the time domain. 

This windowing effect produces a broadening of the peaks 
in the power spectrum and has the added effect of creating a 
greater dynamic range, from the largest spectral component 
to the smallest that can be observed. 

The spectral estimate at the frequency points normally 
is of the form of sin x/x and it is expected that the window¬ 
ing approach will mainly reduce the side lobe effect. When 
the spectral estimate is changed due to the windowing it is 
usually referred to as a modified spectrum in lieu of a raw 
spectrum. 

The need for narrow spectral estimates must be balanced 
with that of the desire to achieve statistical stability. 
Whether or not we choose to use such spectral windows is 

determined by the type of signal we are attempting to analyze. 

2 8 

Brigham has mentioned the idea of noise-like data and 
signal-like data. Statistical stability would play a large 
part in the analysis of noise-like data while practically 
none in signal-type data. 

Modifying the final spectral coefficient Pj^ has the 
effect of smoothing the spectral plot, in the program devel¬ 
oped the smoothing function was an attempt to increase the 
relative smoothness of the plot in order to offer greater 
clarity of presentation so that the sinusoids could be more 
easily distinguished. ; 
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CHAPTER IV 


RESULTS OBTAINED 
Accuracy of the System 

1) Frequency Measurement 

A plot of the system's accuracy as a function of 
the frequency measured in percent of sampling rate is 
shown in figure 26. It can be seen that a limiting 
factor is the size of "DEL", which is the size of the 
spectral estimate. It has been shown that "DEL" is 
directly related to the sampling rate and inversely 
related to the number of data points used. It is well 
to keep in mind that once these parameters are chosen 
the size of "DEL" is fixed. Below a certain level 
of frequency measured compared to the sampling rate, 
the size of "DEL" is a major consideration if our 
error is to be kept in a tolerable range. Above a 
value of 0.1 fs the error was observed to be bounded by 
an upper limit of approximately 0.6%. For a sampling 
rate of 10 khz and 512 data points "DEL" is 19.53 hz. 

It can be seen that if we are attempting to measure a 
500 hz sinusoid we could conceivably incur an error of 
2%. This is seen to increase rapidly as the frequency 
of the sinusoid decreases. For a frequency below 100 hz 
the error is theoretically of considerable size. 

The effect that tends to keep the error fairly con¬ 
stant above the 0.1 fs point is that, as frequency in¬ 
creases, the parameter f/fs tends to approach the Nyquist 
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Rate, 0.5 fs. This decreases the precision of the sam¬ 
ple measurement per cycle, and has been shown to lend 
itself to greater error when the nvimber of scimples por 
cycle approaches a bound of two. Thus the decreasing 
effect of "DEL" is offset by the sampling error. One 
must attempt to choose the sampling rate so that this 
error is minimized in the spectral area of interest. 

2) Amplitude Measurement 

Figure 27 is a plot of the amplitude of the spec¬ 
tral estimates versus peak value of a 1 khz sinusoidal 
input. The digital results show close agreement with 
a standard wave analyzer. It must be noted that, in 
general, the main interest in spectral estimation is 
in relative measurements rather than absolute measure¬ 
ments. Figure 28 is a frequency plot of a 1 khz 
"Half Ramp" function of magnitude 16.5 v. It shows 
good agreement with standard Fourier Series repre¬ 
sentation. 

Signal/Noise Results 

Figures 29 through 32 show the time waveform of a 3.002 ■ 
khz sinusoid of 5.0 volts RMS. White Guassian noise of 0, 

10 v, 20 V 30 V RMS was superimposed on this signal. Figures 
33 through 37 show the frequency spectrum plots of 
these waveforms. 

An Elgenco model 603A noise generator was used as the 
source. This model has a frequency spectrum which is unifbrm 
to plus or minus 0.5 db from 10 to 35000 hz and has a maximum 
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5 Volt Sinusoid 20 Volt Noise 
Figure 31 



5 Volt Sinusoid 30 Volt Noise 
Figure 32 
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spectral density of approximately 2 x 10.5 (volts)^/hz. A 
switch setting was provided to limit the noise bandwidth to 
20 khz. In order to use this noise in testing our system to 
the levels desired it was necessary to use an operational 
amplifier in the Comcor 5000 with a gain of 10. This was due 
to the nominal output range of the generator being 0-3 RMS. 

The noise input was filtered above 4.0 khz in order to avoid 
violating the Nyquist Rate of . 5.0, khz. The concept of White 
Noise is very important./., Wlfi^Nolse has a Guassian Prob- 
ability Distribution and ...a specstral density N{f) that is uni¬ 
form over all frequencies. Noise is considered to be white 
with respect to any given system if the bandwidth of the sys¬ 
tem is small compared to the uniform noise bandwidth. 

Due to the use of the filter and the stated bandwidth 
characteristics of the noise generator our system met these 
requirements. 

For noise with a zero mean value it can be shown that 

= /“ N(f)df = e^ 

2 . • 

where cr = standard deviation 

2 

N(f) = spectral density in volts /hz 
e = RMS value indicated . 

For our system the spectral content of the generator was con- 
stand until 20 khz, therefore 

N(f) = e /20 khz . 

For our amplitude-versus-frequency plot one need only consider 
n(f)^. Therefore as a first assumption: 
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N(f)^ = v/Rf^ 

/2 10 ^ 

= »/rms io“^ 

/2 

The spectral estimate of each point "DEL" is approximately 
10 hz , thus 

N(f)^ = RMS lO”^ volts 

TT 

By increasing the RMS value of the noise in steps of 1 

volt, the net effect should be one of increasing the noise 

value of each estimate by 

k 1 

N(f ) 1 lO"-*- volts = .035 volts 
2 2/2 

Note: The 2 is to account for both positive and negative 
frequencies in the estimate. 

As a conclusion we may estimate that steps of 10 v RMS 
increase the noise level per estimate by approximately 0.1 
volts. Table 2 gives the noise level and signal-to-noise 
ratio attained for the various noise values applied to the 
5-volt RMS sinusoid: 

Table 2 


Noise RMS 

S/N 

Noise Level 

0 V 

29.3db 

0.1120 V 

10 V 

25.8db 

0.1675 V 

20 V 

22.3db 

0.2689 V 

30 V 

19.7db 

0.3460 V 
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The maximum value of the sinusoid in the plot remained con¬ 
stant. The subsequent loss in signal-to-noise ratio is di¬ 
rectly attributable to the increase in noise level and tends 
to generally agree with the theoretical noise calculations 
given above. The average increase per 10 v RMS noise was 
,085 V vice .10 v theoretically. Exact correlation cannot 
be achieved due to the non ideal filter characteristics and 
the finite duration of our sample (.1024 seconds). The noise 
generator is a close approximation to Guassian Noise in the 
practical case only when observed over a time period of suf¬ 
ficiently long duration. Identical measurements were made 
with a HP 302A Wave Analyzer tuned to the frequency of the 
sinusoid. Noise was noted as rapid fluctuations of the 
pointer on the RMS meter. Table 3 tabulates the results 

Table 3 

Voltage Indication (RMS) 

4.80 ± 0 

4.80 ± 0.2 

4.80 ± 0.3 

4.80 ± 0.5 

It should be noted that a noticeable improvement in sig¬ 
nal to noise ratio was achieved in our transition from the 
time domain (-12db) to the frequency domain (19.7db). 
Trade~of£ Considerations 

inherent in the use of any pregram, such as the one ob¬ 
tained, is that the potential user is faced with certain con¬ 
straints which dictate the available courses of action. The 


achieved: 


Noise RMS 
0 

10 

20 

30 
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overriding constraints were limited core storage and the maxi 
mum frequency of sampling. Obviously if one intends to 
tigate all but the low audio range this program is not poten 
tially useful. The maximum sampling rate, as noted earlier, 
is 12.5 khz. This restricts our bandlimitihg condition to 
6.25 khz. The core storage constraint is imposed by the 
equipment capability. In the program developed the total 
core storage available for data was 1025 words. In fact both 
constraining conditions are a direct result of the type of 
equipment used. Having accepted these constraints it is up 
to the user to optimize the amount and type of information 
to be attained, using these constraints as an upper bound. 

The number of samples required is a function of the size 
HARM (Fast Fourier Transform). A minimum number of samples 
must be met, 512 for one version and 1024 for the other ver¬ 
sion. If one is willing to accept the results, as shown 
previously, of quantization to less than a full word, this 
data length minimum may be exceeded many times over. Too 
coarse a quantization scheme, while greatly increasing the 
number of samples taken, may impose an unacceptable bias on 
the system and distort the result, thus voiding any practical 
use. In this way one may achieve a time history of the fre¬ 
quency plot. 

In general one is interested in the precision of measure¬ 
ment available. It has been shown that the larger version of 
HARM offers the greater resolution in sample-point bandwidth 
by increasing the number of points in our sample. This res¬ 
olution is also improved by decreasing our sampling rat^. 
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The greater the number of sample points and the lower the sam¬ 
pling rate the better the resolution, as indicated by the de¬ 
creasing value of DEL . This resolution capability has been 
shown to relate directly to the system's accuracy, thus im¬ 
posing a practical lower limit on the ratio of our particular 
frequency area to the sampling rate. 

If one is to "have his cake and eat it" a fine balance 
must be achieved among these sometimes conflicting paraimeters. 
In essence, the choice is determined by the end results we 
wish to achieve. The final determination must optimize the 
amount of information we wish to obtain by tailoring the 
parameters of the system to the properties of the signals 
under study. 

Ability to Distinguish Between Sinusoids of the Same Amplitude 

The Fourier Transform of a pure tone is in theoiry a 
Delta Fxanction. Since in real life nothing is ideal, a spec¬ 
trum analyzer cannot present an infinitesimally narrow ver¬ 
tical line. Instead the signal is broadened into a pulse 
as shown in Figure 38, a plot of a 3.0 khz sinusoid. Sim¬ 
ilarly, multiple signals close together cause the width of 
the pulses to tend to blend as in Figure 39, which is the 
plot of two sinusoids of equal amplitude at 3.0 khz and 3.1 
khz. These frequencies have been measured by means of a HP 
5245L, digital counter, mentioned previously in this report. 
This illustrates a basic spectrum-analyzer parameter which 
must be considered. The smallest frequency difference be¬ 
tween two equal-amplitude signals that can be displayed is 
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normally defined as the resolution of the analyzer. Two equal 
amplitude signals are considered resolved when they are far 
enough apart to cause a 3 db dip to appear between the two of 
them. 

Larson and Singleton in their report^^, which dealt with 
spectral analysis performed on equipment similar to that used 
in this report, have indicated that for tones of equal ampli¬ 
tude a minimum separation of three resolution cells is needed 
to give a dip between peaks. A resolution cell, as defined 
in their report, is twice the sampling rate divided by the 
number of coefficients determined, that is, 2fg/n, where n 
is the number of data points in the version of HARM used. 

Thus for a sampling rate of 10 khz and the 1024-point version 
of the Fast Fourier Transform a resolution cell would be 9.77 
hz. and the minimum resolution would thus be approximately 
30 hz. 

Figures 39 through 41 confirm their findings and indicate 
additional results as equal-amplitude sinusoids are decreased 
from a 100 hz spacing to 25 hz spacing. It is noted that at 
a 50 hz spacing the threshold effects of merging can be noted. 
The dip is still clearly visible with the level at approx¬ 
imately -15db. At a spacing of only 25 hz, however, the sig¬ 
nals have almost completely merged with the 3.0 khz signal no 
longer clearly distinguishable, its value being only slightly 
greater than the dip between the two tones. 

It may be noted also that one effect has been to increase 
the value of the spectrum coefficient above its previous ap¬ 
proximate value of 5.0. A time waveform is shown for' 
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comparison purposes in Figure 42- For these results the sinu¬ 
soids were set at equal values of nearly 10 volts peak-to-peak. 
The program was no longer capable of distinguishing between 
the two tones and the additional energy from the lower-fre¬ 
quency tone was placed at the 3.025 hz tone with a resultant 
increase in amplitude. It must be noted for clarity that un¬ 
like the Singleton report this is not a power spectrum but 
rather a plot of the discrete Fourier coefficients. In order 
to achieve a power spectrum one must square the amplitude, 
/Cn/^. 



Two Sinusoids of Equal Amplitude Separated by 25 hz 

Figure 42 


It was founa that when the two sinusoids did begin to 
merge it was, as a general rule, the higher-frequency sinu¬ 
soid that tended to predominate. The basis for this result 
in theory has not been fully ascertained but one is prone to 
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speculate that with the presence of noise inherent in any sys 
tern, the sample values achieved would be of a higher value and 
tend more to narrow the apparent frequency difference of the 
two signals. 

A frequency spectrum plot is shown in Figure 43 to demon 
strate that the ability to distinguish between tones is ex¬ 
tendable to n frequencies if the spacing requirement of three 
resolution cells is met. 

General Errors in the System 

Some note must be made of errors in the system that are 
neither inherent in the sampling process nor a function of 
the Fourier Transform subroutine used. The method and type 
of filtering used to bandlimit our signal are very important. 

In theory we would prefer to have an infinitely steep skirted 
filter; however, this is only approached as a goal. The per¬ 
fect filter is assumed to have flat attenuation characteris¬ 
tics over the entire passband. If the concept of the meas¬ 
urement to be performed is one of extreme accuracy, then some 
provision should be made in a supplementary program to take 
into effect the non-ideal filter such as a linear attenuation 
correction over the passband to restore, as nearly as possible, 
the results to the actual case. 

The lack of ideal skirts can cause some "aliasing" if 
we cut our passband too close to the Nyquist Rate. Regard¬ 
less of the value of the cutoff frequency, some distortion 
will result from the attenuation'of the filter skirt and 
should be considered in the final analysis of the results. 
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It will be shown later that failure to use the entire 
capability range of the analog-to-digital converter can se¬ 
verely reduce the dynamic range. It must also be understood 
that when a bit size other than 24 bits is used, the grain 
of our size will be increased and slight variations caused 
by high-frequency components will be attenuated. Thus one 
must xinderstand that for the increased core-storage capabil¬ 
ity some degradation in final presentation must be accepted. 
Dynamic Range 

The dynamic range of a system is normally defined as 
the ratio of power levels, the high level being that of a 
sinusoidal signal at which distortion at some point in the 
system just begins, and the low level being that of the sys¬ 
tem noise with no signal present. It is normally defined in 
terms of decibels, one tenth of the logarithm to the base 
10 of a power ratio: 

Attenuation in decibels = 10 log — 

Most commercially-available Spectrum Analyzers claim a 
dynamic range of 80 db. This however is normally quoted in 
presentation scale capability on the display scope. Larson 
and Singleton have derived a formula for estimating dynamic 
range of a digital—analysis system based on the number of 
bits in the sampling converter and the version of the Fast 
Fourier Transform Subroutine used. Their estimated signal 
to noise rati'o of the spectrum in db is: 

Signal/noise = 3(2K+m) - 1.25 
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where 


K number of bits in the converter 

and 

ni - the value in 2 data values for a spectral 
estimate . 

Thus for the converter used in this system, a 14-bit conver¬ 
ter, (one bit for sign and 13 data bits) and for a 1024-data 
point version of HARM the theoretical signal-to-noise ratio 
was 106.75 db. it is assumed that the full range of the con¬ 
verter is used, in our case 100 volts. 

The maximum signal-to-noise ratio achieved for this sys¬ 
tem was in excess of 60 db. It was achieved on a sinusoid 
of 180-volts peak-to-peak at 1.003 khz as shown in Figure 44. 
The empirical formula for estimation shown also assumes that 
double-precision arithmetic operations are used throughout 
the analysis. For our system this capability was not present. 
The authors stated that the lack of a double-precision capa¬ 
bility would diminish the theoretical signal-to-noise ratio, 
having a much greater effect than lowering the number of bits 
used or the number of data samples in HARM. However, they 
gave no exact relationship. 

The results achieved in this system generally agreed with 
the relationships as set forth in the formula above. An im¬ 
provement was apparent when increasing the number of points 
in the version of HARM. A decrease in the signal-to-noise 
ratio was also noted with a decrease in the number of bits 
used from the converter. These separate effects are easily 
understood in light of previous material covered. Increasing 
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the niJinber of points in HARM has been shown to decrease the 
width of our spectral estimate, thus preventing leakage and 
sharpening the maximim signal levels. Decreasing the nimiber 
of bits used has been shown to increase the grain size of our 
sample and thus increase by a square factor the RMS quantiza¬ 
tion noise. The discrepancy between the basic maximum range 
can be said to be attributable to the lack of double-preci¬ 
sion operations and subsequent roundoff errors that introduce 
system noise. It must also be noted that the amount of 
"windowing" applied to achieve this theoretical relationship 
was not stated. With repeated smoothing one will drasti¬ 
cally reduce the amount of noise level remaining and sub¬ 
sequently achieve a greater signal-to-noise ratio. In the 
results achieved, a smoothing was not applied. The effect 
of the smoothing function can be noted as shown in figure 44 
in comparison to figure 45. Both are plots for the same 
number of data points and converter bits. 

One should attempt to use the full capability of the 
analog-to-digital converter, ± 100 v. It can easily be 
seen that the resolution in amplitude can be increased by a 
factor of 20 db by using the full capability in lieu of a 
10 volt maximum. 
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CHAPTER V 

applicational use op package 

General Conclusions on Use of Package 

The applicational use of digital spectral analysis has 
derived its modern basis from the Fast Fourier Transform. 

The efficiency of such an algorithm has been proven and the 
forming of spectra has been perhaps its most familiar appli¬ 
cation. Specific examples of the analysis of spectra in sig¬ 
nal processing which might be performed digitally with use 
of the Fast Fourier Transform include: seismographic anal¬ 
ysis, speech processing, adaptive digital filtering, gain- 
phase measurements as used in control work with servomech¬ 
anisms, and in general, forming periodograms of random proc¬ 
esses. This is by no means an all-inclusive list but mainly 
is offered to be representative. In recent years a main 
thrust has been in the area of speech processing. ' In 
a guest editorial Bernard Gold aptly stated the challenge, 
"If the speech worker can grasp the theoretical considerations 
behind these sampling and quantizing constraints, he will have 
gone a long way towards reaping the benefits of the digital 
revolution". 

This concept may be applied to all who attempt to use 
such techniques. One must understand the effects of the 
assumption of a "stationary signal", the timing consideration 
imposed by the sampling rate and our bandlimited considera¬ 
tions. Quantizing and subsequent bit-size reduction offer a 
helpful tool in -an attempt to maximize informational content 
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and minimize core storage. They also present a pitfall o 
the unwary if their special auxiliary effects are not noted 
and considered. ' In one application the attendant biasing in 
reducing the bit size may allow a greatly increased sampling 
time, while in other applications it will relegate our results 
to the category of misinformation. The concept of stability 
and bandwidth is of crucial importance when one wishes to make 
the individual spectral estimate as narrow as possible while 
maximizing the estimate's statistical stability. The decision 
to use smoothing is directly related to the type and the qual¬ 
ity of signal we wish to analyze. In general one is faced 
with applying the "trade-off" considerations mentioned pre¬ 
viously to the type of use to which we wish to apply the ' 
system. 

Examples of Applicational Use 
1) Filter Response 

Figure 46 is the frequency spectrum plot of white 
noise passed through a filter. The noise generator was 
an ELGENCO 603A and the filter a Krohn-Hei.tl330A filter 
set at 0.2 hz to 2.0 khz passband. The sampling rate 
was 10 khz 512 data points and a bit size of 12 bits 
was used. The result shown is an average over 8 plots. 
This was achieved by a slight modification of the pro¬ 
gram in order to take advantage of the flatter charac¬ 
teristics of the noise generator if the time of obser¬ 
vation is increased. 

Since this was noise—like data our concern fox 
statistical stability was great. Smoothing was 
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accomplished on the final coefficients. A bit size of 
12 was chosen to increase the core-storage capability 
with negligible distortion observed. The sampling rate 
was chosen at approximately twice the passband in order 
to prevent aliasing from attenuated but sxabstantial 
noise components above the passband of the filter. 

The results shown agree very closely with those 
specified for the Krohn-Hdit filter. The spectrum is 
very flat with a 12-db-per-octave roll-off at the high 
frequency cutoff point. A slight peaking effect is 

1 26 

noted at this point as specified by the filter manual. 

The potential use of this capability of measurement 
in filter design is obvious. The resolution and accu¬ 
racy attained cannot be met with standard spectrum 
analyzers. 

2) Speech Processing 

Much of the present literature in the field of speech 

research today is concerned with speech analysis and 

recognition in terms of a time history of the energy 

versus frequency. Recently comments have appeared in 
32 

the literature critical of the use of Fourier analysis 
in speech recognition. However if one is content to 
analyze the speech waveform with the intent of searching 
for information—bearing elements in the speech process, 
this procedure offers a valuable tool. 
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Figure 47 shows the typical plot of a time history 
of energy versus frequency spectra obtained by analog 
methods. An attempt was made to achieve comparative 
results digitally. 

A male speaker was chosen and three spoken words 
were sampled. The words "Love" and "Hate" were chosen. 
Using "Love" as a reference it is possible to attempt 
to correlate significant aspects of the plot of "Love" 
spoken a second time and "Hate". A bit level of 12 was 
chosen to maximize storage yet avoid undue bias. A 
sampling rate of 10 khz was chosen and a data input of 
512 points to HARM was used. The microphone was input 
into the Krohn-Heit filter, bandlimited from 100 hz to 
3000 hz, the standard voice channel bandwidth. There 
is general agreement in the field that the majority of 
information about speech may be obtained in these ranges. 
The time history plots were separated by approximately 
51 milli-seconds. A 500 milli-second history of the 
monosyllabic words was considered adequate and within 
the limitations of the program. Smoothing was not accom¬ 
plished in these plots, as the author was primarily in¬ 
terested in narrow bandwidths for the spectral estimates 
in the attempt to classify the primary formants in a 
basic recognition scheme. 

A formant frequency is normally defined as that 
frequency which is interpreted fundamentally as the 
frequency of a normal mode of vibration of the human 
vocal mechanism. The term formant frequency, therefore. 
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IS applied both to the normal (or natural frequencies) 
of the vocal system and to their manifestations in the 
speaker s acoustic output, i.e., the frequencies of 
spectral maxima. 

It must be noted .that the analysis was taken in a 
non-pure aco'ustic atmcrsphere; that is, no attempt was 
made to lower the ambient noise level of the IcJsoratory. 
The maximum .amplitude of the voice waveform was limited 
by visual inspection on an oscilloscope at approximately 

V. 

40 volts. 

Figures 48a through 48h show the time history of 
the word "Love" spoken the first time.' Figures 49a 
through 49h show a similar plot for the word "Love" 
spoken a second time. The plots are separated by approx 
imately 50 milli-seconds. The plots of each clearly 
indicate a primary formant for both at approximately 650 
hz with a secondary formant-at approximately 1000 hz. 

One may observe the effect of time on the magnitude of 
each. At approximately 350 milli-seconds the magnitude 
of the formant frequencies have been diminished to an 
extent to be masked by the background noise. The word 
"Hate" in the figures 50a through 50i clearly indicates 
a primary formant at 750 hz, with a secondary formant 
at 2100 hz. At approximately 300 milli-seconds the 
signal has been masked by background noise. 

It is not the purpose of this report to attempt to 
analyze the merits of formant recognition as a means of 
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speech processing and evaluation. However it does show 
the value of this type of analysis as a tool in the 
research environment. 
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CHAPTER VI 


CONCLUSION 

General 

This program for spectral analysis has shown the capa¬ 
bility of a system, based on the Fast Fourier Transform Algo¬ 
rithm, to rapidly and accurately compute reliable power spec¬ 
tra of continuous waveforms. It has shown that a digital 
computer can be effectively used as an analysis tool in the 
frequency domain. It has demonstrated the advantages known 
to be inherent in digital analysis; high resolution, good 
dynamic range and repeatability. The program has shown its 
capability to be used as an adjunct to, rather than a replace¬ 
ment for, analog methods. It is instead a tool to be used 
in a new environment, the digital environment of today's re¬ 
search. It has shown marked advantages over analog methods 
in that the gain and passband stay constant. One result that 
is worthy of consideration is the fact that spectral estimates 
achieved are in digital form and may be immediately used as 
a part in a larger digital problem under investigation. 

Importantly, the program has shown a good general flex¬ 
ibility with a wide variety of possible applications. It 
must also be noted that it served as a valuable tool for the 
author in the development of his capability to use a semi¬ 
hybrid computer facility in an attempt to understand some of 
the theory and application of digital analysis. Moreover, 
it is felt that this package may be a useful tool in further 
research efforts in the audio-frequency area. 
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he results obtained have shown markedly good correla- 

ion with those of similar research being conducted today. 
10,11,14 

Summary of Analytical Results 

The program developed has the capability of sampling at 
rates up to 12.5 khz. Data inputs are available in excess 
of 1024 words. The ability to quantize data results is 
offered in an attempt to allow for an expeinded data capeibil- 
ity so that a time history of the spectra may be developed. 

The frequency-resolution capability has been shown to be the 
sampling rate fs , divided by the number of data points in 
the Fourier Analysis used, n. As an example, for a 10 khz 
sampling rate and 1024 points used, a spectral estimate of 
9.77 hz may be achieved. The accuracy of frequency meas¬ 
urement has been less than 0.2% error in the range of 0.1 to 
0.5 f/fs. Amplitude measurements have shown good agreement 
with standard analog methods and excellent agreement when 
used in relative amplitude measurements of complex wave¬ 
forms. The resolution capability has been shown to be ex¬ 
tremely good, with values of less than 0.3% of the average 
of the signals of interest. Time history spectral plots may 
be achieved separated by less than 50 milli-seconds,'thus 
offering a possible third dimension to the plot, that of time. 
Limitations of the System 

The limitations of the system are those which have been 
imposed by equipment availability and capability. In essence 
the restrictions are those that require a finite time inter¬ 
val and space to achieve the analysis, sampling rate. 
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calculation speed and storage capability. The combination of 
these has limited the analysis range to low frequency (audio 
range) and relatively short duration 

(Time = 1024*24/Nbits*2/fs) 

waveform analysis. The program has been developed so that 
the effect of these limitations on the results achieved may 
be minimized by the ingenuity of the user in optimization of 
trade-offs. 

Recommendations 

The Computer Laboratory is presently in a state of tran¬ 
sition due to equipment modifications and improvement in the 
replacement of the SDS 930 by the SDS 9300. The system capa¬ 
bility has been expanded by the addition of increased core 
storage, display devices, floating-point hardware, and a high- 
capacity drum memory and line printer. 

In general the main recommendation that can be offered 
is that this program be modified to incorporate this system 
into a real-time analysis environment. In order to be real¬ 
istic, the recommendations should only deal with existing 
and proposed equipment modifications. 

The high capacity drum can serve as a temporary memory 
device of access speed far exceeding that of magnetic tape as 
an interface device. The sampling bandwidth will have to be 
reduced in order to do real-time analysis. However with an 
increased storage one may not wish to quantize by varying 
bit size. The author has been able to sample at a rate 
exceeding 16 khz when the need for quantizing has been 


136 



eliminated. The new display device should be used in the 
analysis with an added feature that the output at any time 
may be frozen with a spectral-plot print-out for high 
resolution record purposes. Now that the feasei)ility and 
capability of such a program has been shown^the aspects of 
data smoothing and data overlapping in order to achieve 
closer time separation of plots should bear further inves¬ 
tigation. A closer look at the greater use of the signal¬ 
processing capability of the cinalog computer should be 
investigated and a true hybrid capability approached. 
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$+2 

_166.. 

flRU- 


167 

SKN 

PARAM2 

168 

BRU 

WTING 

-469- - 

-LDA 

iINH2- 

170 

STA 

PARAM2 

171 

BRU 

WTING 

- 172-_ 

— tox. 

svx_ 

173 

LDA 

SV200 

174 

STA 

0200 

- -. 175 

LDA 

SV211 

176 

STA 

0211 

177 

BRR 

INPUT 

_178-CflNT 

--EORM 

.. — 

179 FIX 

form 

9/15 

180 CW 

FIX 

1/CBMM 

181-CfiMM - 

nx 

-UiAZ 

182 

PZE 


183 SPACE 

CBNT 

lSO/0 








184 

TAPIN 

CONT 

ISOiTOTAP 

185 -tec 

RES 

. 1.- . 

186 

ADINT 

BRM 

ADIN 

187 

ININT 

brm 

SAMPLE 

188 

TSTAP 

RES 

130 

189 

TEMP 

RES 

1 

190 

REPEAT 

RES 

1 

_191^MASK-- . 

-QA-TA_OlZZlQQOa 

192 

QNM8K1 

DATA 

040000000 

193 

QNMSK2 

data 

01 

194 

N 

RES 

2 

195 

INPl 

RES 

2 

196 

INP2 

RES 

2 

197 

NUMRCDS 

RES 

2. 

198 

LENGTH 

RES 

2 

199 

SMLENG 

RES 

1 

200 

INP2SM 

RES 

1 

201 

SHIFT 

RES 

1 

202 

SHIFTI 

RES 

1 

203 

SVX - - 


1 





204 

SV211 

RES 

1 

205 

SV200 

RES 

1 

-206 

CeUNT 

RES 

1 

207 

ONTO 

RES 

1 

208 

WDCT 

RES 

1 

-.209 

PARAMl 

-RES. 

1 . . 

210 

PARAM2 

RES 

1 

211 

MA5K1 

RES 

1 

212 

MASK2 

RES 

1 

213 

STORl 

RES 

1025 

214 


END 




AFeRTfiAN US . 

l: SUBROUTINE DATA(FREQ/TIME/NBITS#MljSM8U4Fi/F2) 

2J C THIS SUBROUTINE USES HARM FOR DIGITAL FREQUENCY ANALYSIS 

3J C IT HAS THE CAPABILITY OF CALLING IN DATA FROH THE A/D AND 
4J C analyzing this data AND PLOTTING THE RESULTS 
51 C this subroutine IS LIMITED TO 1024 DATA POINTS 
__FR£Q_1S bur SAMPLING RATE _ - .. 

7: c time is the length of time we wish to sample 
81 c FREQ*TIME must be less than 1024*24/N9ITS 

_ 9: c - Ml specified THE VERSION OF HARM WE WISH -Tfi-USE-- __ 

lo: C 9 being the 512 DATA PT#10 BEING THE 1024 DATA PT VERSION 

11; C SMOU IS 0 FOR NO SMOOTHING USE/I FOR A SMOOTHING FUNCTION 

_121 X _ON OUR FOURRIEILXflEFt . _ 

13: C THIS 6IVED TWO PLOTS,ONE COVERING THE ENTIRE BANDWIDTH 

14: c the other covering only the portion SPECIFIED BY FI AND F2 

151 C ^ FI AND F2 SPECIFY THE DATA BLOW UP WE WANT ON OUR PLOT. 

16: C NBITS IS THE SAMPLE SIZE WE WISH TO USE 

17: DIMENSION A(2048),C(513),F(513),N{128),S(256)/INV(256)/M(3) 

_laj_MCDiMl_ . .... 


19: 

20: 

214_ 

221 

23: 

24: 

M(2)»0 

M(3).0 

_DTiUO/FREQ 

NPTSi»2**M(l) 
NC0 uNT»NPTS/ 128 
DFLii.Q/{DT»NPTS)_ 

25: 

NPL0T»(NPTS/2 )♦! 

26: 

BW«lf0/{DT*2.0) 

271. 

FSTb.0.0 

28; 

XLENG«FREQ*TIME 

29: 

LENGTH«XLENG«ltO 

_aoi_ 

. _NMRCDS»(.XLENGTH*.1)/128)»1 

31: 

NINP1*NBITS*1 

32: 

NINP2*24/NBITS*1 

.334_ 

_NENDiNPLOX.l 

341 

XPLOT»NPLOT 

35; 

EFACTP2.302585 
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IFl^Fl/OEL 

IF2«F2/0EU _ _ - 

IT0T-IF2-IF1 

NREPEAT»(LENGTH+1)*NTIMES/NPTS 
IF(NBITS-14)21,19j19 
19 NRELA-l* 

GO T0 22 

21-N RE LA« N BITS. - __ 

22 THSN»3fO*(2#NRELA*M{l)>-1.25 
IF{M(1).9)38/38/37 

37 NUM.-<LENGTH*l)/1024 
G9 TB 81 

38 NUM»(l.ENGTH*l)/512 

at FJU:WR*4CWt(2.*»23)-li ) .. ___ 

PRINT28/FREQ/BW/DEL#NPTS/TIME 
PRINT30/NBITS 

PRlNT25a/NINPl/NlNP2/LENGTH/NMRC0S/NBiTS 
PRINT263/NC0UNT 
FtD-OiO 

-O0- t OO- U 2/NPLftT- -...— 

561 F(I)*FST+DEL 

.571 FST»F(I) 

- 5ai— lOO-CBNttNUE -- _ _ 

59) FSTnO.O 

60) REWINDl 

-att- .- ..CALW- tr4RU T (Ntfa p i/NINP2/LENGTH/ NARCOS/NBfm— - 

62) c this is the meta-symbbl callable FBRTRAN subrbutine 

63) REWINDl 

. _ hht- kk Oft-3 I*1/NC0UNT 

65) READ TAPEl/{N(J)/J*l/128) 

66) DS 3 vi-1/128 

-674-tlxt2»j*t)--*lt-U*256- - 

68) IIP1>IU1 

69) A{II)«FLBATF(N(J))*FACTBR 

-7W-A.< MPU »0^0- - 

71) 3 CBNTINUE 

72) CALL HARM(A/M/INV/S/-1/IFERR) i 


36 ) 

-^74 

38) 

39) 

40) 

41) 

42) 

- 434 - 

44) 

45) 

- -464 - 

47) 

48) 

- _43: 

50) 

51) 

- 52! 

53) 

54) 
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73J 

7U___ 

751 

76: 

77: 

78: 

79: 

aoj_ 

8i: 

82: 

S3: 

84: 

851 

.a6J.__ 

87: 

88 : 

89: 

90: 
9i: 


99: 

loo: 

_iQir- 

102: 

103: 

105: 

106: 

108: 

109: 


DO 10 iPl^NPTS 

_ 

I2M1pI2.1 / 

A(1)»SQRT((A(I2)**2)+a(I2M1)»*2)/ 
10 CONTINUE 

1F(SM0U*1.0)59^49/49 
49 PRINT265 
-50 D0.51^Jj2/NENDl — . 

KP1pI*1 

KM1pI*1 

Cl 1 )f (•A(KMl}44f0tA{n«A(KPl) }/4tQ 
IF(C(I)-0.0)101/51/51 
101 C(I)«OtO 

__51-CfiNIlNUE_.... 

GO TO 70 
59 PRINT264 
6Q-D0-61/li2/NEND , 

61 C(I)»A(I) 

70 C(1)»A(1) 

_ClNELmiAINEUaiJ__ _ . _ 

CAV«0*0 

CT0T»0.0 

_CMAXaQ*0 _ - „ 

DO 18/I*1/NPL0T 
IF(CMAX«C(I) )16/16/13 


13 .CT0T»CT0T*C(I) 

18 CONTINUE 
. _^PRINT266/CMAX . 
CAViCTOT/XPLOT 
DO 78/KS1/2 

_CREEJtQjiD-- - . ... 

DO 73/Iil/NPUeT 
IF(C(I)*CAV)72/72/71 
7Jl_CREPi£REP*C( I) - 
GO TO 73 

72 CREP^CREP+CAV 





noi 73 CONTINUE 

-liil. -- 70 CAVaCREP/XPtOT .... . 

U2l ST0N«20.0*AL0G(CMAX/CAV)/EFACT 

U3| PRlNT2S$4STeN 

llAJ PRlNT260iTHSN 

1151 PRlNT987#CAV 

116S CALL F4PL0T(F>C*NPL0T) 

-14.71 _ . na 33i».i«uiTaT __ _ 

1181 iriN«lTlFl 

119J A(I)«C(IFIN) 

120 J F4U«F4IFIN) 

1211 333 continue 

122; IFISENSE SWITCH3)11/12 

123V-C — T HlS■ S EN S&-SWITCH IS USED FOR DETERM1>HJJG--WHETHER OR .NB-T- 

124! C TO PRINT OUT THE FOURRIER COEF. OF INTEREST 
125! 11 PRINT57 

126: PRINT27# ( (Fd )/A(l > ) / I»l/ ITOT) 

127! 12 CALL F4PL0T(F/A/ITOT) 

128! 27 F0RMAT(4(lX/0PFlOt2/2X/lPlEl0.3)) 

._i2ai--28-Ffl RM A T . aHl * 1 0X/ tSAMPLING RATEit/F-lOta»10X4lBWt«rlO>X/.5Xx»Q£L4i8-^ . 

130! •F10.2/5X/*NUM PTS • $/F5.0 #$TIME»$/ F8»4 /IX/tSECONDS*) 

131! 30 FORMAT!11X6USING A BIT LEVEL OF $/ 15) 

4421 40 FORKATUIX^PLOT «/I4/* OF t/IA) . 

133! 57 FORMAT{//$FREOUENCY AND FOURRIER COEF* •) 

134! 258 FORMATdlX/ tlNPlit/I3/2X/$INP2^$/14, 2X/fLENGTHf •/18, 2X 

44S4-•NBITS*4^;4)- ---- 

136! 259 FORMAT!11X/4SIGNAL TO NOISE RATIO OF PLOT • t/FStl/f DBtl 

137! 260 F0RMAT!11X/*THE0RETICAL S/N OF PLOT * /F5»l) 

138! 263 FORMAT!11X/»THE NUMBER OF RECORDS FOR ONE INPUT TO HARM4/I10) 

1391 264 FORMAT! llX/tUSING NO SMOOTHING FUNCTION$.) 

140! 265 FORMAT!llXtUSlNG A SMOOTHING FUNCTION*) 

IkU _266 .FftfiMAXd4X*3MAX. VALUEii^ 1P1E10.3) ___ . 

142! 987 F0RMAT!11X/4N01SE LEVEL IS */lPlE10.3) 

1431 999 FORMAT!llX/fPLOT */I5/*0F •/I5) 

_l4iAl- r eturn- - --- 

1451 END 

••GLOBAL variables*^* 






APPENDIX B 


SDS 

930 COMPUTER MACHINE 

INSTRUCTION LIST 

CODE 

MNEMONIC 

NAME 

00 

HLT 

HALT 

01 

BRU 

BRANCH UNCONDITIONALLY 

02 

EOM 

ENERGIZE OUTPUT M 

06 

EOD 

ENERGIZE OUTPUT TO DI¬ 
RECT ACCESS CHANNEL 

12 

MIW 

MEMORY INTO W BUFFER 
WHEN EMPTY 

13 

POT 

PARALLEL OUTPUT 

14 

ETR 

EXTRACT 

16 

MRG 

MERGE 

17 

EOR 

EXCLUSIVE OR 

20 

NOP 

NO OPERATION 

23 

EXU 

EXECUTE 

32 

WIM 

W BUFFER INTO MEMORY 
WHEN FULL 

33 

PIN 

PARALLEL INPUT 

35 

STA 

STORE A 

36 

STB 

STORE B 

37 

STX 

STORE X 

40 

SKS 

SKIP IF SIGNAL NOT SET 

41 

BRX 

INCREMENT INDEX AND 
BRANCH 

43 

BRM 

MARK PLACE AND BRANCH 

46 

;- 

REGISTER CHANGE 

50 

SKE 

SKIP IF A EQUALS M 

51 

BRR 

RETURN BRANCH 
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CODE 

MNEMONIC 

NAME 

52 

SKB 

SKIP IF M AND B DO NOT 
COMPARE ONES 

53 

SKN 

SKIP IF M NEGATIVE 

54 

SUB 

SUBTRACT 

55 

ADD 

ADD M TO A 

56 

sue 

SUBTRACT WITH CARRY 

57 

ADC 

ADD WITH CARRY 

60 

SKR 

REDUCE M, SKIP IF NEG 

61 

MIN 

MEMORY INCREMENT 

62 

XMA 

EXCHANGE M AND A 

63 

ADM 

ADD A TO M 

64 

MUL 

MULTIPLY 

65 

DIV 

DIVIDE 

66 

— 

SHIFT RIGHT 

67 

— 

SHIFT LEFT 

70 

SKM 

SKIP IF A = M ON B MASK 

71 

LDX 

LOAD INDEX 

72 

SKA 

SKIP IF M AND A DO NOT 
COMPARE ONES 

73 

SKG 

SKIP IF A GREATER THAN M 

74 

SKD 

DIFFERENCE EXPONENTS 

AND SKIP 

75 

LDB 

LOAD B 

76 

LDA 

LOAD A 

77 

EAX 

COPY EFFECTIVE ADDRESS 
INTO INDEX REGISTER 
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APPENDIX C 


BASIC ALGORITHM USED 
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Given three integers M^.» = 1, 2, 3 and a three- 
dimensional array of complex numbers Afki, ko. ko) 
where k^ = 0 , 1 .n^,-! and , this sub¬ 

routine computes the three-dimensional, complex 
Fourier series: 


N,-l N^-l N 3 -I 


XOi. jg. Jg)"* 2 2 2 A(k ,k ,k )e 

iTlo rto 12 3 


where 


2t1I 


^2’h h-h 


N, 


N„ 




( 1 ) 


j = 0,1.Ny-1 and - 1,2,3. 


The inverse series may also be computed: 



Ng -1 

Ng-l 

A(ki, kg, kg) - N ♦N *N 2 

2 

2 

ii=o 

J2=0 

^3=° 


X(kj^, k^, kg) e 




N, 


N„ 


N„ 


For one dimension, the algorithm for the complex 
Fourier transform of complex numbers A (k) is: 


( 2 ) 


N-1 

Xa)=Z^ A(k)w^‘' (3) 

k= 0 

where w= e-~—, N= 2^, j-= 0, 1, 2, ..., N-1 

expresses X (j)as a fimction of the M arguments 
-1, j M - 2 .. Jl' io binary representa¬ 

tion of j: 


. •• oM-1 


+ jj2 + 


Jn. J. 


0 or 1 
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1 1 ' 


xa 


M 


. • • • • • IJ, IJ^) ™ ^ ^ ^ 1.^ 

VO 1^1=0 1^1-1= 

A(kM-l.V'^V 


,, 41r 0»-l 

Since J “Si -1 - v/ **Sl -and 1. 

the iimermost Bum In (4) yields an array: 


^1^0*'Si-2.V 


E 

■St-i' 




M-l 


Successive arrays .,Aj^ obtained by summing 

over kj ^_2 .k^ respectively are computed. The 

general formula is: 

.^L-I’Vl-I.V 


“ S “^L-l ^0.^L-2’ 'Sl-L. 

'Sl-L=’ 0 




( 6 ) 


The final array will bo the desired X, The storage 
Indexing convention used is to let the M arguments 
of AlCJq, ..,., kg) be the binary representation of 
the index of the storage location for AL(jo,... ,ko). 
In this way, each step of the algorithm involves 
fetching from two storage locations and returning 
results in the same two locations, thereby saving 
storage. However, the elements of the final array: 

.Vl' 

are in the wrong order. The program reverses the 
order of bits in bit representations to compensate. 
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